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ABSTRACT 
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(N ■ 

We investigate the global transition from a turbulent state of superfluid vorticity 
(quasi-isotropic vortex tangle) to a laminar state (rectilinear vortex array), and 
vice versa, in the outer core of a neutron star. By solving numerically the hydro- 
dynamic Hall-Vinen-Bekarevich-Khalatnikov equations for a rotating superfluid 
>- ■ in a differentially rotating spherical shell, we find that the meridional counterflow 

driven by Ekman pumping exceeds the Donnelly-Glaberson threshold through- 
out most of the outer core, exciting unstable Kelvin waves which disrupt the 
rectilinear vortex array, creating a vortex tangle. In the turbulent state, the 
torque exerted on the crust oscillates, and the crust-core coupling is weaker than 
in the laminar state. This leads to a new scenario for the rotational glitches ob- 
o served in radio pulsars: a vortex tangle is sustained in the differentially rotating 

q ■ outer core by the meridional counterflow, a sudden spin-up event (triggered by 

an unknown process) brings the crust and core into corotation, the vortex tangle 
relaxes back to a rectilinear vortex array (in 10 5 s), then the crust spins down 
electromagnetically until enough meridional counterflow builds up (after £ 1 yr) 
^ ■ to reform a vortex tangle. The turbulent-laminar transition can occur uniformly 

or in patches; the associated time-scales are estimated from vortex filament the- 
ory. We calculate numerically the global structure of the flow with and without 
an inviscid superfluid component, for Hall-Vinen (laminar) and Gorter-Mellink 
(turbulent) forms of the mutual friction. We also calculate the post-glitch evo- 
lution of the angular velocity of the crust and its time derivative, and compare 
the results with radio pulse timing data, predicting a correlation between glitch 
activity and Reynolds number. Terrestrial laboratory experiments are proposed 
to test some of these ideas. 
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1. Introduction 

Timing irregularities in a rotation-powered pulsar, such as discontinuous glitches (Lyne 
et al. 2000; Zou et al. 2004; Shabanova 2005) and stochastic timing noise (Hobbs 2002; Scott 
et al. 2003), provide an indirect probe of the internal structure of the star. The physical 
processes usually invoked to explain these phenomena are (un)pinning of Feynman-Onsager 
vortices in the crystalline inner crust (Anderson & Itoh 1975), starquakes (Ruderman 1976), 
and thermally driven vortex creep (Alpar et al. 1984b; Link et al. 1993). Less attention 
has been directed at the global hydrodynamics of the superfluid, except within the context 
of the spin-up problem in cylindrical geometry (Anderson et al. 1978; Reisenegger 1993; 
Carter et al. 2000). The importance of the global hydrodynamics was demonstrated by 
Tsakadze & Tsakadze (1980), who simulated pulsar rotational irregularities in the laboratory 
by impulsively accelerating rotating containers of He II, obtaining qualitative agreement with 
radio timing data (e.g. glitch amplitudes and post-glitch relaxation times). 

In this paper, we examine how the global flow pattern of superfluid in the outer core of a 
neutron star affects the rotation of the star. We focus on the outer core for simplicity: vortex 
pinning is thought to be weak or non-existent there (Alpar et al. 1984a; Donati & Pizzochero 
2003, 2004), and the fluid is mainly isotropic ( 1 S'o Cooper pairing) (Sedrakian & Sedrakian 
1995; Yakovlev et al. 1999), reducing the problem to a hydrodynamic one in a spherical shell. 
Even with this simplification, the calculation remains numerically challenging: the spherical 
Couette problem for a superfluid was solved for the first time only recently (Peralta et al. 
2005), generalizing previous work on the cylindrical Taylor-Couette problem for a superfluid 
(Henderson et al. 1995; Henderson & Barenghi 2004) and the spherical Couette problem for 
a classical viscous fluid (Marcus & Tuckerman 1987a; Dumas & Leonard 1994). 

An isotropic ^So-paired) neutron superfluid is described by the two-fluid Hall-Vinen- 
Bekarevich-Khalatnikov (HVBK) model. In this model, the viscous normal fluid and inviscid 
superfluid components feel a mutual friction force whose magnitude and direction depends 
on the distribution of Feynman-Onsager vortices (Hall & Vinen 1956a; Bekarevich & Khalat- 
nikov 1961). The vortices are organized in a rectilinear array if the flow is strictly toroidal, 
but they evolve into a tangle of reconnecting loops when the counterflow along the rota- 
tion axis exceeds a threshold, exciting the Donnelly-Glaberson instability (DGI) (Glaberson 
et al. 1974; Swanson et al. 1983; Tsubota et al. 2003). Peralta et al. (2005) showed that the 
DGI is excited in a neutron star under a wide range of conditions, driven by the meridional 
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component of the spherical Couette flow (SCF) in the interior. The mutual friction force 
changes dramatically during transitions between a vortex array and a vortex tangle (Gorter 
& Mellink 1949; Vinen 1957; Swanson et al. 1983; Schwarz 1985), affecting the rotational 
evolution of the star. 

In this paper, we propose a phenomenological model for timing irregularities in radio 
pulsars based on the creation and destruction of a vortex tangle — superfluid turbulence 
(Barenghi et al. 1995; Vinen 2003) - - in the outer core of a rotating neutron star. In 
our scenario, a glitch comprises the following sequence of events, (i) Differential rotation 
between the outer core and crust of the star, built up over time through electromagnetic 
spin down, generates a meridional Ekman counterflow in the outer core. We show that the 
axial counterflow exceeds the DGI threshold, creating a vortex tangle throughout the outer 
core. The mutual friction in this turbulent state takes the isotropic Gorter-Mellink (GM) 
form and is much weaker than the mutual friction associated with a rectilinear vortex array, 
(ii) When the glitch occurs, triggered by an unknown mechanism, the outer core and inner 
crust suddenly come into corotation and the vortex tangle decays, ultimately converting into 
a rectilinear vortex array. The decay process lasts ^ 10 5 s, depending on the drag force 
acting on the vortex rings, after which the mutual friction takes the anisotropic Hall- Vinen 
(HV) form (Hall & Vinen 1956a,b) and increases by ~ 5 orders of magnitude, precipitating 
a "torque crisis", (iii) After the glitch, differential rotation builds up again between the 
outer core and the crust due to electromagnetic spin down. When the axial counterflow 
exceeds the DGI threshold, after ^ 1 yr, the vortex array breaks up again into a tangle and 
the mutual friction drops sharply. Similar transitions from turbulent to laminar flow in a 
superfluid have been observed in laboratory experiments where He II, cooled to a few mK, 
flows around an oscillating microsphere (Niemetz et al. 2002; Schoepe 2004). 

The paper is organized as follows. HVBK theory is briefly reviewed in §2, together 
with the pseudospectral numerical method which we use to solve the HVBK equations in a 
differentially rotating spherical shell. The physics of the turbulent-laminar transition in a 
generic glitch scenario is elaborated in §3. The response of the stellar crust to a turbulent- 
laminar transition in the outer core is calculated numerically in §4. Finally, the results are 
summarized and applied to observational data in §5; a fuller observational analysis will be 
carried out in a future paper. 

2. HVBK model of the outer core 

The neutron superfluid in the outer core exists in an isotropic ( 1 5'o) phase for densities p 
in the range 0.6 < p/p* < 1.0 (where p* is the nuclear saturation density) and an anisotropic 
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( 3 P 2 ) phase for 1.0 < p/p* < 1.6 (Yakovlev et al. 1999); the depth at which the %- 3 P 2 
transition occurs is not known precisely (Epstein 1988). As the hydrodynamic equations for 
the 3 P 2 phase (Mastrano & Melatos 2005) are complicated and hard to handle numerically, 
we construct our model around the 1 Sq phase, described by HVBK theory in this paper. We 
plan to extend the model to the 3 P 2 phase in the future. 

The protons in the outer core are probably in a type II superconducting state, where the 
magnetic field is quantized into fluxoids (Sauls 1989). Protons inside magnetic fluxoids in- 
teract with neutrons in Feynman-Onsager vortices (Sauls 1989; Ruderman 1991; Link 2003). 
If vortex pinning occurs, and if the core magnetic field has comparable poloidal and toroidal 
components (Thompson & Duncan 1993), the geometry of the neutron-fluxoid interaction 
can be complicated (Ruderman et al. 1998). To avoid these difficulties, we incorporate the 
protons (and other charged species) into the normal component of the superfluid, a common 
approximation (Comer & Joynt 2003; Prix et al. 2004). We neglect proton-neutron entrain- 
ment (Mendell 1991; Sedrakian & Sedrakian 1995), which is observed in terrestrial 3 He- 4 He 
mixtures (Andreev & Bashkin 1976; Andersson & Comer 2001) and is an important mecha- 
nism modifying the mutual friction for rectilinear vortices in a neutron star (Mendell 1991; 
Andersson et al. 2006). We neglect hydrodynamic forces arising from entropy gradients and 
angular momentum textures (Mastrano & Melatos 2005). Finally, we neglect vortex pinning 
for simplicity, except at the inner and outer boundaries, even though it seems likely that 
pinning plays an important role in glitch dynamics (Epstein & Baym 1988; Baym et al. 1992; 
Hirasawa & Shibazaki 2001; Link & Cutler 2002). 1 



2.1. HVBK theory 

The HVBK model for a rotating superfluid is a generalization of the Landau-Tisza two- 
fluid model that includes the hydrodynamic forces exerted by quantized vortices in He II 
(Hall & Vinen 1956a; Bekarevich & Khalatnikov 1961) and 1 S , -paired neutron matter (Tilley 
& Tilley 1986; Prix et al. 2004). Fluid particles in the theory are assumed to be threaded 
by many coaligned vortices, a valid assumption over length-scales longer than the average 
inter- vortex separation (and shorter than the average radius of curvature). In the continuum 
limit, the vorticity of the superfluid component satisfies w s = V x v s ^ macroscopically 
(cf. V x v s = microscopically). The normal component comprises charged species (protons 
and electrons, locked together by the external magnetic field) plus thermal excitations and 



1 Three-fluid models analyzing post-glitch relaxation favor pinning at the core boundary (Sedrakian & 
Sedrakian 1995), e.g. due to vortex cluster-Meissner current interactions (Sedrakian & Cordes 1999). 
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behaves like a classical, viscous, Navier-Stokes fluid (kinematic viscosity v n ). The isothermal, 
incompressible (V • v n = V • v s = 0) HVBK equations of motion take the form (Barenghi & 
Jones 1988; Henderson & Barenghi 2000) 

—7— = hz/ n Vv n H F Vw a , (1) 

at p p p 

^ = -Y2± + v,T-^F-™V\ U .\, (2) 
at p p p 

with d n:S /dt = d/dt + v„ iS • V, where v n (v s ) and p n (p s ) are the normal fluid (superfluid) 
velocities and densities respectively. Effective pressures p s and p„ are defined by Vp s = 
Vp — \prNiyns) an d Vp„ = Vp + |p s V(v^ s ), with v ns = v n — v s . Note that, in neutron star 
applications, the gravitational potential can be subsumed into p without loss of generality, 
in the incompressible limit. 

The vortex tension force per unit mass, which arises from local circulation around quan- 
tized vortices (Andronikashvili & Mamaladze 1966), is defined as 

u s T = w s x(Vx u s ), (3) 

with uj s = u) s /\uj s \. Here, v s = (k/Att) ln(fe /oo) is the stiffness parameter, k = h/2m n is the 
quantum of circulation, m n is the mass of the neutron, do is the radius of the vortex core, and 
60 is the intervortex spacing. Note that u s , which has the dimensions of a kinematic viscosity, 
controls the oscillation frequency of Kelvin waves excited on vortex lines (Henderson et al. 
1995). 

The mutual friction force per unit mass, F, arises from the interaction between the 
quantized vortex lines and the normal fluid (via roton scattering in He II and electron 
scattering in a neutron star) (Hall & Vinen 1956a,b). The structure of this force depends on 
the global configuration of the vortices. If the vortices form a rectilinear array, the friction 
takes the anisotropic HV form (Hall & Vinen 1956a,b) 

F = l -Bu s x (u> s x v ns - u s T) + \b'{u s x v ns - v s T) (4) 

where B and B' are temperature-dependent mutual friction coefficients (Barenghi et al. 
1983). If the vortices form a vortex tangle, the friction takes the isotropic GM form (Gorter 
& Mellink 1949) 

— '(^)v„., (5) 
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where A' = 5 3 p^7r 2 x?/3p 2 X2 is a dimensionless temperature-dependent coefficient, related 
to the original GM constant (usually denoted by A in the literature) by A' = Apu. Here, 
Xi and X2 are dimensionless constants of order unity (Vinen 1957). The HV and GM forces 
are in the ratio ~ 10 5 under the physical conditions prevailing in the outer core of a neutron 
star (Peralta et al. 2005). This implies that the normal and superfluid components of the 
star (and hence the crust and core, through viscous torques) are effectively uncoupled most 
of the time (when enough differential rotation has built up), but become tightly coupled in 
the immediate aftermath of a glitch when a turbulent-to-laminar transition occurs (Peralta 
et al. 2005). It is important to bear this in mind when reading §3. 



2.2. Spherical Couette flow 

We model the global hydrodynamics of the outer core by considering an HVBK super- 
fluid enclosed in a differentially rotating spherical shell, i.e. superfluid spherical Couette 
flow. The inner (radius Ri) and outer (radius R 2 ) surfaces of the shell rotate at angular 
velocities Q± and f2 2 , respectively, about a common axis. All quantities are expressed in 
dimensionless form using R 2 as the unit of length and fi^ 1 as the unit of time. We adopt 
spherical polar coordinates (r, 9, (f>), with the z direction along the axis of rotation. 

The boundary conditions satisfied by the normal and superfluid components are sub- 
tle. Neither component can penetrate the boundaries, implying (v n ) r = (v s ) r = 0. The 
tangential normal fluid velocity satisfies the no-slip Dirichlet condition, viz. (v n )# = 0, 
(y n )<t> = Ri, 2^1,2 sin 9. The behavior of the superfluid at the boundaries is influenced by the 
quantized vortices, which can either slide past, or pin to, irregularities on the boundaries. 
The former scenario implies (uJ s )e = (w s )^ = 0, while the latter implies (<jj s ) r <jj s x v L = 
and requires a knowledge of the vortex line velocity Vl locally. As the HVBK model pro- 
vides no information about v^, we adopt the widely used no-slip compromise (v s )g = 0, 
( v s)<?i — ^1,2^1,2 sin 0, which in the spherical case does not restrict the orientation of the 
vortex lines at the surface (Henderson et al. 1995; Henderson & Barenghi 1995, 2004). Note 
that, initially, v n and v s must be divergence-free, with V x v s ^ 0; the Stokes solution 
with v s = v n satisfies this constraint (Landau & Lifshitz 1959). However, the numerical 
results presented in §3 and §4 are obtained after several Ekman times (and hence rotation 
periods) have elapsed, by which time the flow has forgotten its initial conditions, as in most 
astrophysical applications. 

Many existing glitch models are based on cylindrical geometries (Anderson et al. 1978; 
Alpar et al. 1984b; Abney & Epstein 1996; Larson & Link 2002), except for a few experi- 
ments with spheres (Tsakadze & Tsakadze 1973, 1980). However, the global flow pattern in 
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spherical Couette flow differs from its cylindrical counterpart in two important ways. First, 
in a cylinder, the principal flow is toroidal, whereas, in a sphere, the principal flow is an 
axisymmetric combination of toroidal flow and meridional circulation for all Re (Tuckerman 
1983). 2 Second, it is misleading to approximate the equatorial region of a sphere by a cylin- 
der, even for fli « i?2; the curvature, although slight, causes significant differences in the 
critical Taylor number at which vortices appear (Soward & Jones 1983; Stuart 1986). The 
differences between cylindrical and spherical Couette flow are equally prominent in classical 
fluids and superfluids (Peralta et al. 2005). 

Spherical Couette flow of a viscous fluid is controlled by three parameters: the Reynolds 
number Re = VL\R\jv n , the dimensionless gap width S = (R?,— R\)/ R2, and the angular shear 
AQ — fi 2 — Jli- If a superfluid is also introduced, v s emerges as an additional parameter. In 
the laminar regime (Re £ 10 4 ), the global flow can be classified according to the number of 
cells of meridional circulation in each hemisphere. In a narrow gap (5 < 0.11), the merid- 
ional circulation is slow and can be approximated by (v n ) e /(R 2 AQ) 10~ 2 5 2 (i? 2 /-Ri)-fte 
(Yavorskaya et al. 1986). As Re increases, a centrifugal instability generates toroidal Taylor 
vortices near the equator until, at the onset of turbulence (Re 10 4 ), a helical traveling wave 
develops. In a wide gap (5 > 0.4), by contrast, all secondary flows are nonaxisymmetric, 
and they oscillate when Re exceeds a threshold (Yavorskaya et al. 1986). In classical fluids, 
these flows have been thoroughly investigated numerically and experimentally (Yavorskaya 
et al. 1975, 1977; Belayev et al. 1978; Yavorskaya et al. 1986; Marcus & Tuckerman 1987a,b; 
Junk k Egbers 2000). 

Recently, high-resolution numerical simulations of superfluid spherical Couette flow were 
performed successfully for the first time (Peralta et al. 2005). 3 Figures la and lb display 
streamlines of the normal fluid and superfluid components for Re = 3 x 10 4 , 5 = 0.3, 
Afl = 0.1, and HV mutual friction. Both components exhibit a similar structure in each 
hemisphere: a "square" Ekman cell near the pole, three to four secondary meridional cells, 
a number of smaller cells, and two equatorial vortices near the outer boundary which are 
created and destroyed intermittently. We emphasize that the normal component is not 
assumed to be in uniform (solid body) rotation; it is affected by the superfluid (through 
mutual friction) and the boundary conditions (differential rotation). Figures lc and Id 
display streamlines for a run with identical parameters and GM mutual friction, for which 



2 Spherical Couette flow is a combination of quasi-parallcl-plate Ekman pumping at the pole and quasi- 
cylindrical, ccntrifugally driven Taylor vortices at the equator (Junk & Egbers 2000). 

3 Visualizing superfluid spherical Couette flow experimentally is challenging. It has been attempted suc- 
cessfully using glass microspheres in turbulent He II (Bielert & Stamm 1993). The technique of particle 
image velocimetry is currently being evaluated (Zhang et al. 2004). 
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the inter-fluid coupling is ~ 10 5 times smaller, as noted in §2.1. We observe that the normal 
fluid exhibits a similar number of vortical structures, which are better developed and more 
regular than in Figure la, and the superfluid flow pattern is nearly unchanged from Figure 
lb. Meridional circulation is apparent in all the results; we find axial velocity components 
in the range 1(T 5 <, (v n ) 2 , (v s ) z £ 1(T 2 . 

The flow pattern in Figure 1 changes with temperature just like in cylindrical Couette 
flow. The superfluid behaves like a classical, viscous fluid near the critical temperature. The 
Taylor vortices elongate axially, with a more complex pattern of eddies and counter-eddies 
emerging, as the temperature is lowered (Henderson et al. 1995; Henderson & Barenghi 
1995). For Re £j 268, anomalous modes emerge in the normal component, as in classical 
axisymmetric flow (Lorenzen & Mullin 1985), characterized by Ekman cells rotating in the 
opposite sense to the classical flow (Henderson & Barenghi 2000) except near the critical 
temperature. For Re £ 268, the tension becomes less important than the mutual friction and 
the superfluid component comes to resemble the normal component (Henderson & Barenghi 
1995). 

3. Glitch-induced turbulent-laminar transition 

In this section, we investigate how the vorticity state of the outer core changes before, 
during, and after a glitch, in the context of a generic glitch scenario where the trigger 
mechanism of the glitch is unspecified. 

3.1. Before a glitch: vortex tangle 

Laboratory experiments on the attenuation of second sound in narrow channels (Vinen 
1957; Swanson et al. 1983), and numerical simulations based on the vortex filament method 
(Tsubota et al. 2003), show that, in a rotating container of He II, an axial (along z) coun- 
terflow v ns > t>DG = 2(2f2 2 z / s) 1 / 2 excites growing Kelvin waves which destabilize a rectilinear 
vortex array (Glaberson et al. 1974). The Kelvin waves grow exponentially until reconnec- 
tion between adjacent vortices occurs, generating a dense tangle. 4 For typical neutron star 



4 The Kelvin waves excited by the DGI are unrelated to the Kelvin waves generated by oscillations of the 
nuclear lattice in the inner crust (Epstein & Baym 1992). 
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parameters, one finds 

n n 1/2 



taking ln(fe /oo) = 20. 



UDG = L56 I To^k^ ' cms_1 ' (6) 



The axial component of the Ekman counterflow in the outer core of a typical neutron 
star generally exceeds the instability threshold (6). This is illustrated by Figure 2, a greyscale 
plot of (v ns ) 2 /i>DG; the DGI is active in regions with (v ns ) z /t> DG > 1. Figure 2a corresponds 
to HV mutual friction (initially rectilinear vortex array); Figure 2b corresponds to GM 
mutual friction (after a tangle forms). In both figures, the DGI is active in most of the 
computational domain, implying that an initially rectilinear array is disrupted and, once 
disrupted, stays that way. This result is extracted empirically from the simulations. The 
inclusion of compressibility and hence stratification (Abney & Epstein 1996) restricts the 
DGI region to a thin boundary layer, but the overall conclusion is the same (see §4.5). 

The axial counterflow that excites the DGI is driven by Ekman pumping. As the crust 
spins down electromagnetically at a rate Q*, differential rotation builds up between the crust 
and outer core, with 

Asi - = 3 16 x 10-6 ( io-'L-' ) fe) rad '" 1 ' (7) 

where t is the time elapsed since the last glitch (Lyne & Graham-Smith 2006; Lyne et al. 
2000). The differential rotation induces meridional circulation (Greenspan 1968; Reisenegger 
1993): an Ekman boundary layer, with (v n )# ~ R*A£l, develops on a time-scale ~ 27r/fi*, 
and grows radially to a thickness d E ~ Re~ x l 2 R* cm on a time scale t E ~ iJe 1 / 2 ^*) -1 , 
spinning up the interior fluid (Andersson 2003; Andersson et al. 2005). Only the normal 
fluid is Ekman pumped directly; the superfluid is spun up by the normal fluid through the 
mutual friction force (Adams et al. 1985; Reisenegger 1993). 

In order to estimate the counterflow velocity, we compute (v ns ) 2 empirically from nu- 
merical experiments and scale up the results to neutron star parameters. For example, in 
a typical set of runs with 5 = 0.3 (HV mutual friction), 0.4 (HV mutual friction), and 0.5 
(GM mutual friction), and Re = 3 x 10 4 , we find (v n ) 2 ~ 3.15(v s ) z , (v n ) 2 ~ 0.78(v s ) 2 and 
(v„) 2 w 0.06(v s ) 2 , respectively, on average over the grid, implying (v ns ) z ~ (v n ) 2 ~ R*AQ cm , 
as a general rule and hence 
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In a typical neutron star, the Reynolds number in the outer core (temperature T) is 
very large, viz. 



given the standard viscosity v n resulting from electron-electron scattering (Flowers & Itoh 
1979; Cutler & Lindblom 1987; Andersson et al. 2005). Consequently, the flow is likely to be 
turbulent (Alpar 1978). Our numerical simulations cannot access this regime due to com- 
putational limitations (Peralta et al. 2005). However, we know from laboratory experiments 
with He II that the flow at Re ^ 10 5 closely resembles classical, Navier-Stokes turbulence: 
the superfluid and normal fluid vorticity evolve in concert, because the vortex lines are locked 
to the normal fluid eddies by mutual friction (Ashton & Northby 1975; Barenghi et al. 1997). 



For the DGI to be triggered in superfluid turbulence, two conditions must be met. 
First, there must be vortex segments directed parallel to the counterflow, i.e. uj s ■ w ns ^ 
(Barenghi et al. 1997). As the superfluid is locked to the normal fluid, and the vortex density 
is high, the condition u) s ■ v ns ^ is satisfied locally throughout the core. Second, a vortex 
segment is locally unstable to (helical) Kelvin wave perturbations of wavelength A if the 
local counterflow exceeds the critical speed t>jxj b = k ln(A/27rao)/ (2A) (Barenghi et al. 1997). 
Large wavelengths (of the order of the size of the star) are most unstable, with A ~ 10 6 
cm giving v^ G b = 3 x 10~ 2 cm s -1 . In this regard, the DGI threshold (6) is conservative 



During a glitch event, the differential rotation instantaneously shuts off and the tangle 
starts to decay. In a self-sustaining vortex tangle, Kelvin waves grow until their amplitude 



5 In He II experiments, it is possible to reach a vortex density ~ 10 7 cm 2 (Skrbck et al. 2000), comparable 
to that in a neutron star, viz. 10 5 (f2*/10 2 rads -1 ) cm -2 (Baym et al. 1992; Lyne & Graham-Smith 2006). 

6 A vacillating flow (Barenghi et al. 2004) can also trigger vortex oscillations by inducing an oscillating 
counterflow. A perturbation with wavelength A grows if the vacillation frequency oj is mantained over a 
time-scale t < \v ns /2iri' s uj (Barenghi et al. 2004). 



Re 




(9) 



5 



(v$g<V DG ). 6 



3.2. Glitch: decay of the tangle 
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approaches the average vortex separation, whereupon the vortex lines continuously reconnect 
to form loops (Jou & Mongiovi 2004). Reconnection stops once the counterflow ceases, as 
happens immediately after a glitch, when the outer core and inner crust come into corotation 
(Qi = f^)- The decay time-scale equals, to a good approximation, the reconnection time- 
scale just before the counterflow ceases. 

To estimate r^, we approximate the vortex loops by rings of radius R, whose character- 
istic lifetime is given by = R 2 (2u s a)^ 1 , with a = Bp n /2p (Barenghi et al. 1983; Tsubota 
et al. 2004). The radius of a ring can be estimated from the vortex line density L (length 
per unit volume), with R = 0.5L -1 / 2 in a steady-state tangle. The evolution of L is given by 
Vinen's equation (Vinen 1957), generalized by Jou & Mongiovi (2004) to include rotation: 



1 dX 
th~dt 



-a 3 X 2 + 



OllV r , 



x 3/2 



(3 4 v r , 



X, 



■" (10) 

with X = kL/Q 2 , where a iy Pi are dimensionless friction coefficients, whose values depend 
on the temperature (Mongiovi & Jou 2005). By fitting to experimental data, one finds 
«3 = 20.0ai, /3i = 35.6ai, fa — 53.6«i, and /?4 = 1.43ai (Swanson et al. 1983; Jou & 
Mongiovi 2004). From experiments, we also have a = Xi a i an d Xi ~ 1 (Vinen 1957), so we 
take a = ai below. Setting dX/dt = in (10) and solving for L in the steady state, we find 
that the tangle decays on a time-scale 



2.5 x 10 3 /t 
af„ s ln(6 /ao) 



7.6 x 10 5 ( 



a 



10- 7 



-i 



-2 N _ 2 

\ ft 



X I lO-^rads-'J ylyrj B ' (12) 

where v ns is the counterflow speed immediately before the differential rotation shuts off, 
given by (8). 

Estimates for the friction parameter a are uncertain (Sonin 1987), but calculations based 
on electron scattering by proton vortex clusters and magnetized vortex cores give B 10~ 4 
(Sauls et al. 1982; Sedrakian & Cordes 1998). 7 This implies a < 10~ 7 for p/p n = 10~ 2 in 



7 An absolute lower limit is B = 10 18 , calculated from electron scattering off vortex clusters in the high- 



the outer core (Sedrakian & Sedrakian 1995). Therefore Td is greater than the glitch trigger 
time-scale, as assumed a priori. Note that the B and B' parameters were also calculated by 
Mendell (1991) and Andersson et al. (2006), who found B' = B 2 and 



B 



4 x 10 



x 




7/6 



(13) 



where x p = pj p v is the proton fraction, and m p and m* are the bare and effective masses of 
the proton respectively (Andersson et al. 2006). Equation (13) includes modifications from 
entrainment (Mendell 1991; Andersson et al. 2006). In the outer core, with p = 2.8 x 10 14 
g cm -3 and x p = 0.038 (Mendell 1991), we obtain B = 2.7 x 10 -4 , approximately the value 
quoted by Sauls et al. (1982) and Sedrakian & Cordes (1998). 

In laboratory and numerical experiments, it is observed that the vortex tangle is polar- 
ized; the average vorticity projected along the rotation axis is not zero (Finne et al. 2003; 
Tsubota et al. 2004; Tsubota & Kasamatsu 2005). This is because, under these experimen- 
tal conditions, one has Pi ^ P^v ns j '(k^) 1 ^ 2 in equation (10), where Pi is the polarization- 
inducing term. Naively, one might expect the tangle to be even more polarized in a rapidly 
rotating neutron star, but in fact the opposite is true: the tangle is less polarized, because 
we have Pi <C PkV ns j \kVL 2 ) 1 ^ 2 from Figure 2. A thorough study of this issue, including 
whether the large-^2 limit considered by Jou & Mongiovi (2004) in deriving equation (10) 
is applicable for v ns ^> (kQz) 1 / 2 , lies outside the scope of this paper. 



After the tangle decays to a rectilinear array, differential rotation due to electromagnetic 
spin down between the inner and outer shells begins to grow until (v ns ) z exceeds foe- From 
equations (6) and (8), this occurs after a time 
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density regime, where the relaxation time-scale equals the age of the pulsar [Sedrakian & Cordes (1998); A. 
Sedrakian 2006, private communication]. 



whereupon a vortex tangle develops again in the outer core, via the DGI. This triggers a 
switch from HV to GM mutual friction, weakening the coupling between the normal fluid 
and superfluid components. 

The characteristic time-scale for the tangle to develop equals the growth rate of helical 
vortex perturbations (Kelvin waves) via the DGI. By linearizing Schwarz's equation (Schwarz 
1985, 1988), in the vortex filament model, one finds that the fastest growth occurs at a 
wavenumber k = v ns /2v s , with (Tsubota et al. 2004) 

/cln(& /a ) 

~ 4.9xio 4 f^-rY A r Vs, 

\10~ 7 J \10 2 rads- 1 J 

where (16) follows from (15) by substituting (6). Hence, for typical neutron star parameters, 
the tangle is reestablished over ~ 1 d after ~ 1 yr elapses following a glitch, less than the 
inter-glitch interval observed in most pulsars (Shemar & Lyne 1996; Lyne et al. 2000). Note 
that r g does not equal t^, the time required for a tangle to grow from a rectilinear array is 
shorter than the time required for an existing tangle to decay back to a rectilinear array, 
provided that t < 8.0t ta n (and longer otherwise). Note also that v ns just before the glitch 
(~ R*AQ), which appears in (11), is typically greater than v ns at t tan (~ ^dg)) which appears 
in (15), provided that t > t tmi . 

Equation (16) gives the minimum time for a tangle to reform, assuming it does so 
simultaneously everywhere in the outer core of the star. In practice, v ns does not exceed 
v-dg everywhere simultaneously. Regions where the DGI is activated are interspersed with 
regions where it is not, even at r = Ri and R2, so that the transition from HV to GM friction 
is more gradual than equation (16) suggests. This important issue, which is also relevant to 
the interpretation of the laboratory experiments performed by Tsakadze & Tsakadze (1980), 
is explored briefly in §5. However, a thorough numerical study lies outside the scope of this 
paper. 



(15) 
(16) 



4. Rotational evolution of the crust 

In this section, we calculate numerically the torque acting on the stellar crust during 
a transition from a state of turbulent vorticity (tangle) to laminar vorticity (rectilinear 
array), by solving the hydrodynamic HVBK equations in a differentially rotating spherical 
shell. Global simulations of this type are necessary to calculate the observable response of a 
neutron star, e.g. Q(t) and Cl(t), to the internal physics elaborated in §3. 
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4.1. Numerical method and parameters 

We solve equations (1) and (2) using a pseudospectral collocation method to discretize 
the problem (Boyd 2001; Canuto et al. 1988) and a time-split algorithm to step forward 
in time (Canuto et al. 1988), as described by Peralta et al. (2005). The velocity fields are 
expanded in a restricted Fourier series in 6 and and a Chebyshev series in r (Orszag 
1974; Boyd 2001), with (JV r , N , iV ) = (120, 250, 4) modes required to fully resolve the flow. 
The initial time-step At = 10~ 4 is lowered to At = 10~ 5 during the GM — > HV transition 
to prevent spurious oscillations in the torque (Peralta et al. 2005). Instabilities arising 
from the sensitivity of spectral methods to the boundary conditions (Peralta et al. 2005), 
and oscillations due to the Gibbs phenomenon (Gottlieb et al. 1984; Canuto et al. 1988), 
are smoothed using a low-pass spectral filter (Don 1994), a common practice (Osher 1984; 
Mittal 1999). 

We adopt parameters as close to those of a realistic neutron star as our computational 
resources permit. In the outer core, we take p s /p = 0.99 and p n /p = 0.01 (Baym et al. 
1992), and hence B = 1.5, B' = 0.9, and A' = 1.0 x 10~ 4 at the corresponding scaled 
temperature T/T c in He II (Barenghi et al. 1983; Donnelly & Barenghi 1998). These friction 
coefficients are ~ 10 4 times greater than those used in the analytic estimates in Section 
3. We adopt the higher values deliberately. Otherwise, the effects introduced by B and 
B' would take too long to build up; computational limitations prevent us from simulating 
more than ~ 10 rotation periods. (For the same reason, our AQ/Q is unrealistically large.) 
Realistic Reynolds numbers, estimated from equation (9), are too high to be simulated 
directly, so we restrict ourselves to the range 10 3 ^ Re ^ 10 5 ; the most stable evolution 
occurs for Re ~ 3 x 10 4 , corresponding to an Ekman number (Re)' 1 ~ 3 x 10 -5 that is two 
orders of magnitude smaller than in a typical neutron star (Abney & Epstein 1996). The 
tension force is dominated by mutual friction, except in relatively old (> 10 4 yr) neutron 
stars (Greenstein 1970), with v s ~ 10~ 18 in our dimensionless units. For 10 < A < 10 4 , we 
find empirically that the torque depends weakly on the GM force, so we use 10 < A < 10 2 
to obtain identical results at lower computational cost. 



4.2. Turbulent-laminar transition 

We consider a turbulent initial state with Qi > Q 2 and (v ns ) z > vbg everywhere, 
such that the normal and superfluid components are coupled by GM mutual friction. The 
transition from a vortex tangle to a rectilinear vortex array is simulated by changing the 
mutual friction suddenly, from GM to HV, and simultaneously spinning up the outer shell, 
such that f2 2 = This occurs at t = 20 in all the figures to follow. The angular shear 



- 15 - 



before the transition lies in the range 0.1 < AQ < 0.3, for gaps in the range 0.2 < 5 < 0.5. 
Note that we are forced to choose Afl unrealistically large, compared to typical neutron star 
values, in order to allow the inner and outer surfaces of the shell to "lap" each other within 
a reasonable run time. 8 Various combinations of AQ and 5 affect the viscous torque on the 
inner and outer surfaces during steady differential rotation in different ways, as discussed 
below. 

Streamline snapshots for a run with 5 = 0.5 and AQ = 0.2 are presented in Figures 
3a-3f, with the turbulent-laminar transition GM — > HV, ^2 —> ^i occurring at t = 20. The 
HV mutual friction, with |Fhv|/|Fgm| ~ 10 5 , couples the normal and superffuid components 
strongly. 9 Both components exhibit a quasi-periodic motion, with a secondary circulation 
cell filling most of the shell, and two or three smaller vortices emerging intermittently. An 
oscillatory polar Ekman cell also exists initially, disappearing at t = 140. 

We can compare this behaviour to the flow without a turbulent-laminar transition. 
Figures 4a-4f display a sequence of streamline snapshots for 5 = 0.5 and Ail = 0.2, with 
GM mutual friction at all times. At t = 18, three polar cells are observed in both components, 
together with a secondary cell at mid-latitudes. The latter cell elongates after t = 20 and is 
present only in the superfluid component at later times. At the equator, the flow switches 
between one and three cells, but this behavior persists after t = 100 only in the superfluid 
component, while the normal fluid decouples increasingly. 

The axial counterflow is significant at all times in Figures 4 and 5. Scaling to neutron 
star parameters, we find v ns /vr>G ~ 10 3 (Af2/10~ 4 rads~ 1 )(r2*/10 2 rads~ 1 ) near the equator, 
and v ns /v DG ~ 10 6 (Afi/10" 4 rads- 1 )(fi*/10 2 rads- 1 ) near the poles, at t = 200. 

Figures 5a-5d display AQ/Q as a function of time after the turbulent-laminar transition 
for 5 = 0.2, 0.3, 0.4, and 0.5, respectively. 10 The evolution is qualitatively similar in all 
the cases considered. There is an initial transient, in which the flow adjusts to the initial 
spin up on a time-scale t ~ 6, just like in the classical (viscous) spin- up problem (Greenspan 
1968). The fractional torque deviation, Ati/Cl, is plotted as a function of time in Figures 
6a-6d. An initial jump, coinciding with the sudden spin up, is observed in all the panels, 



8 Flow structures associated with differential rotation develop on the time-scale (Ail) 1 . 

9 The empirical result |FhvI/|Fgm| ~ 10 5 is supported by order-of-magnitude estimates (Peralta et al. 
2005) but should be treated with caution, because the microphysics of the GM friction force (and hence the 
value of A') has not yet been worked out rigorously in a neutron superfluid. 

10 Note that Afi/fi* is not calculated self-consistently: it is obtained by numerically integrating the viscous 
torque on the outer sphere, after solving (1) and (2) subject to the boundary condition f2 2 (*) = constant, 
i.e. we do not include the (slight) evolution of f^W when computing the fluid flow. 
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followed by a rapid exponential decay and persistent, small- amplitude oscillations. We find 
an e -1 relaxation time of r ~ 2, 1.6, and 1.6 for AQ = 0.1, 0.2, and 0.3 respectively, with 
r essentially independent of 5. In contrast, in the inverse experiment, where the mutual 
friction changes from HV to GM, Peralta et al. (2005) observed that Ati is almost constant 
after the post-glitch transient. The period of the oscillations does not change with t for 
t > 100. 



4.3. Effect of the superfluid 

The effect of the superfluid on the evolution can be analyzed with the help of Figure 
7a, which plots the torque on the outer sphere as a function of time (solid curve) in the 
absence of spin up (5 = 0.5, Re = 3 x 10 4 , AQ = 0.1, GM friction). For comparison, 
the dashed curve plots the torque for a classical viscous fluid with the same parameters 
mentioned above. The evolution is qualitatively similar in both cases, displaying oscillations 
with the same periodicity and similar amplitude, but the curves diverge for t ^ 20 as the 
viscous fluid evolves more rapidly to a stationary state, driven by classical Ekman pumping. 
Classical theory predicts a time-scale t E ~ 173 in the limit where the Rossby number AQ/Q 
vanishes; an exponential fit to the dashed curve gives ~ 243, which agrees well given that 
we have AQ/Q = 0.1 in our numerical experiments. In a superfluid, on the other hand, 
the torque decays faster initially, with time constant ~ 50, then starts to increase again for 
t ~ 240, as part of a long-period oscillation whose ultimate fate is unknown as it occurs over 
a time-scale beyond our computational limit. 

Figure 7b shows the evolution of the fractional change in angular velocity after spin 
up (f2 2 — 0.9 — > 1.0 at t = 20) in three scenarios: instantaneous change in mutual friction, 
GM — > HV (solid curve); constant GM friction (dashed curve); and a classical viscous fluid 
(dashed dotted curve). The classical viscous fluid and constant GM friction evolve similarly. 
The torque decays exponentially on a time-scale ~ 2.15, then decays almost linearly with t. 
By contrast, in the GM — * HV transition, the torque decays exponentially, on a time-scale 
~ 2.20, then oscillates with peak-to-peak amplitude w 0.4 (in units of pR\Vl\) and period 
~ 6. In other words, the oscillations are sustained by the HV friction, which is much weaker 
than its GM counterpart. Note that the oscillations do not correspond to Tkachenko waves, 
since the HVBK equations assume that the free energy of the vortex array depends only on 
the vortex line density, not on vortex lattice deformations (Chandler & Baym 1983, 1986; 
Donnelly 1991). 
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4.4. Streamlines 

The meridional streamlines for the transitions described in the previous paragraph evolve 
as in Figure 8. Figures 8a-8b show the streamlines for a classic viscous fluid before (t = 20) 
and after (t = 22) the glitch. There is a rapid redistribution of vorticity near the outer 
sphere, with two circulation cells replaced by one elongated meridional shell, while the flow 
near the inner sphere hardly changes, with one circulation cell near the equator and one 
near the pole. Figures 8c-8d show the meridional streamlines for the normal component 
of the superfluid. The flow resembles a classical viscous fluid (cf. Figure 8a), although it 
does become more complex after the spin up. As shown in Figure 8d, the large cell near the 
outer sphere is replaced by one primary circulation cell and three secondary cells (including 
one near the poles). Note that the greatest contributions to the torque come from the mid- 
latitude regions where the structure of eddies and counter-eddies is richer, as in Figure 8d. 
Finally, in Figures 8e-8f, we see the meridional streamlines for the inviscid component of the 
superfluid. The pattern before the jump (Figure 8e) differs from the normal fluid; the GM 
friction creates additional eddies near the outer sphere and closer to the poles. However, 
after the transition to HV friction, the inviscid component comes to resemble the normal 
component, as can be seen by comparing Figure 8f with Figure 8d. This occurs because the 
two components are coupled more strongly by HV friction (|Fhv|/|Fgm| ~ 10 5 ), so that the 
superfluid is dragged along by the normal fluid. 

4.5. Stratification 

Gravitational stratification can strongly suppress Ekman pumping, as the Ekman layer 
is squashed close to the outer sphere (Clark et al. 1971; Abney & Epstein 1996). Nevertheless, 
the Ekman layer, however thin, always exists, in order to satisfy the boundary conditions 
at r = R 2 . It contains a meridional counterflow given by equation (8). Consequently, 
the transition from turbulent to laminar vorticity (and vice versa), which relies on such a 
meridional counterflow, still occurs in a stratified star, and its dramatic effect on the torque 
is the same. The Ekman layer, no matter how thin, acts like a film of "oil" between two 
sliding surfaces (here, the outer core and inner crust), whose coupling strength ("stickiness") 
changes abruptly when the meridional speed of the "oil" exceeds a threshold. In other words, 
equations (8), (12), (14), and (16), and the conclusions that follow from them, are unaltered 
by stratification; even though the volume of the outer core occupied by the vortex tangle 
(which does not affect the torque on the crust directly) does change. 

The Ekman layer thickness decreases as e~ KY d E , where d E is the thickness without 
stratification and Ky is the compressibility of the fluid. The Ekman time-scale increases as 
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(e Ky — l)/ty 1 tE, where £e is the time-scale without stratification. The fluid is restricted to 
move on concentric spherical shells (Levin & D'Angelo 2004). Importantly, (e Ky — l)^ 1 ^ 
is the time for the Ekman layer to extend throughout the outer core, not the time required 
to establish the meridional flow at r as R 2 (~ ^i 1 ), which controls the onset of the DGI. 

It is conceivable that the DGI is excited in shells at certain radii where v ns peaks, so that 
we end up with a sequence of shells containing alternating laminar and turbulent superfluid 
vorticity. In this scenario, the detailed study of which lies outside the scope of this paper, the 
stability of the vorticity configuration is restricted by the Richardson criterion, which states 
that the configuration is Kelvin- Helmholtz unstable for N 2 \dv,p/dr\~ 2 < 1/4 (Mastrano & 
Melatos 2005). In a neutron star, the Brunt- Vaisala frequency is N ~ 5x 10 2 s" 1 (Reisenegger 
& Goldreich 1992). 

The effects of stratification are not considered in detail in this paper because they cannot 
be studied properly with our numerical method; the pressure projection step only works with 
divergence-free velocity fields (Bagchi & Balachandar 2002; Peralta et al. 2005). However, 
a crude approach to get a feel for the effects is to numerically suppress v r using a low-pass 
exponential filter (Don 1994), viz. v r — > exp[— (k/N r )' y In e]v r , with < \k\ < N r , where 
e = 2.2 x 10 -16 is the machine zero, and 7 is the order of the filter. If we ramp up 7 with r 
as 7 = r(n 2 — ni)/(i? 2 — Ri) + ri\ — R\{n2 — ni)/(i?2 ~ -Ri), with n\ — 2 and n 2 = 12, then 
v r is weakly suppressed near R2 (7 = 12) and strongly suppressed near R\ (7 = 2), so that 
the flow is confined approximately to concentric shells. Preliminary results, for a viscous 
Navier-Stokes fluid, show that the filtering reduces the meridional circulation, flattening the 
streamlines radially. 



5. Summary and discussion 

In this paper, we investigate how transitions between turbulent and laminar states of su- 
perfluid vorticity alter the standard theoretical picture of pulsar rotational irregularities like 
glitches and timing noise, (i) Most of the time, except in the immediate aftermath of a glitch, 
differential rotation in the outer core drives a nonzero, poloidal counterflow which continu- 
ously excites the DGI. A vortex tangle is thereby maintained in the outer core. The mutual 
friction in this regime, which is of GM form, couples the normal and superfluid components 
loosely, (ii) Immediately after a glitch, the differential rotation ceases, as does the poloidal 
counterflow. The vortex tangle decays over the mean life-time of its constituent vortex rings, 
Td = 7.6xl0 5 (a/10~ 7 ) _1 (f2*/10 _13 rads~ 2 )^ 2 (t/l yr)~ 2 s. A rectilinear vortex array develops, 
and the mutual friction switches to HV form, coupling the normal and superfluid compo- 
nents much more strongly, (iii) After t tan = 0.49 (^/lO 2 rads- 1 ) 1 / 2 (^/lQ- 13 rads" 2 )" 1 
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yr, electromagnetic spin down builds up the differential rotation sufficiently to drive a 
poloidal counterffow that exceeds the DGI threshold. A vortex tangle forms again in a 
time r g = 4.9 x 10 4 (a/10 _7 ) _1 (£X/10 -13 rads -2 ) -1 s, and the mutual friction reverts to GM 
form. Note that vortex pinning provides the boundary conditions for the superfiuid SCF 
but does not occur within the outer core itself (Donati & Pizzochero 2003). Therefore our 
new phenomenological picture is not a complete model for glitches. It merely clarifies the 
vorticity state of the outer core before and after a glitch as an input into future models than 
incorporate the full glitch dynamics, including trigger mechanisms related to pinning in the 
inner crust. 

We draw together the strands of the model in Figure 9, which displays the evolution of 
the torque and the regions where the DGI is active during the following numerical experiment: 
we fix Af2 = 0.1 until t = 20, accelerate the outer sphere instantaneously to corotation at 
t = 20, then decelerate the outer sphere according to Q 2 (t) = 1 — 0.001 (t — 20) for t > 20. 
This mimics the situation in a real pulsar, where we have *C t tan , i.e. Ekman pumping 
brings the fluid into corotation before the DGI gradually starts being reexcited throughout 
the outer core. To make the experiment as realistic as possible, we do not assume that the 
mutual friction takes the same form everywhere in the outer core, but rather choose GM 
or HV friction at each point according to whether (v ns ) 2 is greater or less than t> DG locally. 
In this comparison, we approximate v ns by v n , as in equation (8), because v s can become 
very complicated (e.g. Figure 1), creating numerical difficulties. In order to satisfy £e <C £tan 
while keeping Af2 large enough so that the spheres "lap" each other several times, we are 
forced by computational exigencies to adopt a relatively low Reynolds number Re = 100, 
shortening the Ekman time (£e ~ lOfi^ 1 <C ^tan), and to artificially boost i>dg, so that 
|(v n ) z /i>DG| does not exceed ~ 3 throughout the computational domain. 

The results of the above numerical experiment are presented in Figure 9. Contours of 
|( v ri) 2 /^DG|) before and after the spin up at t = 20, are plotted in Figures 9a-9e. Shaded 
regions indicate where the DGI is active, i.e. |(v n ) z /wDG| > 1- Just before the glitch 
(Figure 9a), 32 % of the superfiuid is in a turbulent state, with the DGI active close to the 
inner sphere and at intermediate latitudes where meridional circulation is significant. After 
the differential rotation shuts off at t = 20, the DGI initially spreads through the shell as 
transient axial flows increase, occupying 39 % of the volume at t — 22. However, the flow 
quickly settles down to a state of near-corotation during the time interval 24 £ t £ 50, HV 
friction dominates, and the torque decays exponentially, with time constant ~ 10 (followed 
by a linear decay). At t — 50, the DGI slowly begins to assert itself again, starting from the 
inner sphere. As for a classical viscous fluid, when the outer sphere spins down, it pumps 
fluid radially inward and along the axis of rotation (Vanyo 1993), so the axial flow speed is 
greatest near the inner sphere. By t = 120, when AQ = 0.1, the vorticity state is similar to 
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that at t — 20, before spin up. 

One might wonder if, in a realistic neutron star, the superfluid ever exits the turbulent 
state and becomes laminar. For the simulations in this paper, which have Re < 3 x 10 4 , 
the answer is yes. Figure 9 shows that, when the glitch occurs and the spheres come into 
corotation, |(v n ) z | falls below the DGI threshold after a time t ~ t E - The vortex tangle is 
then guaranteed to decay on a time-scale given by (11) and (12), as observed in terrestrial 
experiments. However, for more realistic neutron star Reynolds numbers (Re > 10 8 ), which 
are too challenging to simulate at present, the turbulent eddies in the normal fluid decay 
more slowly when the spheres come into corotation, so that |(v n ) 2 | remains above the DGI 
threshold for longer. If this happens, the vortex tangle may persist until the next glitch 
occurs, so that the superfluid never exits the turbulent state. 

What are the implications of Figure 9 and the results in §3 and §4 for observations 
of glitches? Before considering this question, we emphasize again that the results in this 
paper do not constitute a theory of glitches, because important questions regarding the 
glitch trigger remain unresolved. Nevertheless, some general remarks can be made. First 
of all, it is clear that transitions between flow (and vorticity) states in superfluid SCF are 
caused by changes in Re, and that such transitions become more frequent and complicated 
as Re increases (Yavorskaya et al. 1977; Junk & Egbers 2000). This is compatible with the 
observation that adolescent pulsars (~ 10 4 yr old, like Vela) glitch most actively (Lyne et al. 
2000). In younger pulsars (age £ 10 4 yr), T and hence v n are relatively high, so Re is low. 
In older pulsars (age £ 10 4 yr), Q and hence Re are low following electromagnetic spin 
down [although this trend is not straightforward and can be masked by localized heating 
from differential rotation between the superfluid and the crust (Greenstein 1975; Larson & 
Link 1999) or crust cracking (Link et al. 1998; Franco et al. 2000)]. A systematic statistical 
study of glitch activity versus Re will be published elsewhere (Melatos et al. 2006), but 
preliminary estimates of T and hence Re from cooling curves (Tsuruta 1974, 1998; Page 
et al. 2004) including superfluidity (Flowers & Itoh 1976, 1979; Andersson et al. 2005) give 
Reynolds numbers in the range 10 s ^ Re ^ 10 12 for glitching pulsars. Two of the most active 
glitchers, the Crab and Vela, have Re ~ 10 9 and Re ~ 10 10 respectively. One expects that, 
at such high Re, the fluid is turbulent, with the kinetic energy concentrated at large scales 
(Yavorskaya et al. 1978, 1986), as for a classical viscous fluid (Smith et al. 1993; Barenghi 
et al. 1997). This suggests that superfluid turbulence in pulsar interiors is an important 
factor in glitch dynamics. 

If it is true that the vorticity in the outer core exists in a turbulent state before a glitch, 
as postulated in our model, then t tan represents a lower bound on the time between glitches. 
In testing whether this bound is respected by the glitching pulsars currently known, we are 
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hampered by the fact that most of these objects have only glitched once. Nevertheless, for 
all the 28 pulsars that have glitched repeteadly, we find that the minimum inter-glitch time 
interval t min is greater than t tan , as the theory predicts (Melatos et al. 2006). The object 
PSR 2116 + 1414 approaches the bound most closely, with t min = 3.9 yr and t tan = 2.0 yr. 
This is encouraging, because the 28 objects cover five decades in t tan and three decades in 
t m in, and the theoretical expression (14) for t tan contains zero free parameters. Note that 
ttan is proportional to the characteristic age (= f^/^iX) divided by Vt l J 2 . Note also that the 
activity parameter defined by McKenna & Lyne (1990) involves glitch amplitudes (which are 
highly variable) as well as mean recurrence times, so we do not predict a correlation between 
the activity parameter and t tan . 

It is harder to test the theoretical decay time-scale of the vortex tangle, as predicted 
by (12), because it remains unclear what observable features are engendered by the decay 
process. The observed exponential post-glitch relaxation is of viscous origin and occurs on a 
time-scale much larger than t&. On the other hand, the decay of the tangle is accompanied 
by a large increase in mutual friction (GM — > HV), which may be connected with the rapid 
jump in Q during a glitch. The jump in Q has never been resolved in time, in pulsars 
which are nearly constantly monitored, consistent with the predictions of (12) for the Crab 
(r d = 3 x 10~ 4 s) and Vela (r d = 0.2 s). However, equation (12) predicts that it may be 
possible to resolve the Q jump in older pulsars, provided that the time between glitches 
does not increase faster than fi*. In making these estimates we assume the canonical value 
a = 10~ 7 for every object, in the absence of a microscopic theory, yet this is clearly an 
oversimplification because a is sensitively temperature dependent. 

Oscillations in (l* were observed before (period ~ 10 d) and after (period ~ 25 d) 
the Vela Christmas glitch, with Ati/Cl* « 0.17 (McCulloch et al. 1990). In our numerical 
simulations, persistent torque oscillations are always observed when the outer core rotates 
differentially, as occurs before a glitch. They are also observed after a switch from GM to HV 
friction, as occurs after a glitch. By comparing the dashed and solid curves in Figure 7b, we 
see that the oscillations are sustained by HV mutual friction. The oscillation period in our 
simulations is much shorter than in pulsar data, because we are restricted to Re < 3 x 10 4 . 
An alternative explanation is that vortices in the inner crust oscillate relative to the normal 
fluid in the core (Sedrakian et al. 1995). 

Several of the effects explored in this paper have been studied in terrestrial laboratories. 
Our results will motivate new experiments of this sort, cf. Alpar (1978) and Anderson et al. 
(1978). Although it is hard to access the neutron star regime p n <C p s in He II, where inter- 
atomic forces are appreciable, transitions between turbulent and laminar superfluid vorticiy 
have been observed in experiments with microspheres immersed in 4 He at mK temperatures 
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(Niemetz et al. 2002). Promising results on the relaxation of rotating He II were obtained 
by Tsakadze & Tsakadze (1980), but again these results are for p n £ p s and hollow spheres 
rather than a differentially rotating shell. We propose to extend these experiments in two 
directions: (i) by investigating low-Rossby-number (AQ/Q <C 1), high- Reynolds-number 
(Re 3> 10 5 ) SCF with He II at the temperature which minimizes p n /p s ; and (ii) by repeating 
(i) with a nonideal dilute-gas Bose-Einstein condensate confined in a differentially rotating 
magneto-optical trap, in order to probe the stability of a vortex lattice to Kelvin wave ex- 
citations in the regime p n <C p s (Parker & Adams 2005). The presence of a vortex tangle 
in He II can be detected by standard second-sound absorption techniques (Hall & Vinen 
1956a; Swanson et al. 1983), and the torque in experiment (i) can be monitored to look for 
oscillations when a change from HV to GM friction (or vice versa) is triggered by the DGI. 

In classical Navier-Stokes fluids, injection of vorticity into a metastable laminar state 
can trigger turbulence, e.g. seed vortices injected into a cylindrical vessel containing 3 He-B 
(with T < 0.6T C ) generate a vortex tangle that eventually decays into a rectilinear vortex 
array (Finne et al. 2003). Unlike He II, the normal component in 3 He-B is laminar in these 
experiments and does not participate in the turbulent dynamics, due to its comparatively 
high viscosity \y n ~ 1 cm 2 s _1 ^> u s , cf. v n ~ v s in He II; see Finne et al. (2004)]. Standard 
glitch theories assume that the normal component is tightly coupled to the crust by the 
external magnetic field (Alpar & Sauls 1988; Jahan Miri 1998). This suggest a second 
possible SCF experiment with 3 He-B (at T < 0.6T C ) in a hollow spherical container in 
which the normal fluid corotates with the container, and therefore does not participate in 
the DGI, while the superfluid is free to become turbulent. Such an experiment can probe 
what aspects of the turbulent-laminar transition are caused by the normal and superfluid 
components respectively. Transitions between a vortex tangle and a rectilinear array can be 
detected using non-invasive nuclear magnetic resonance techniques (Finne et al. 2003). 
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Fig. 1. — Streamlines in superfluid spherical Couette flow with 5 = 0.4, Re = 3 x 10 4 , and 
AQ = 0.1. (a) Normal fluid with GM mutual friction, (b) Superfluid with GM mutual 
friction, (c) Normal fluid with HV mutual friction, (d) Superfluid with HV mutual friction. 
The streamlines are calculated by integrating the in-plane components of the velocity fields 
in the plane x — at t — 18. 

Fig. 2. — Normalized counterflow velocity {^ ns )z/ v DG m superfluid spherical Couette flow 
with 5 = 0.4, Re = 3x 10 4 , and Aft = 0.1 at t = 46. (a) HV mutual friction, (b) GM 
mutual friction. 

Fig. 3. — Meridional streamlines of the normal fluid (left) and superfluid (right) components, 
for 5 = 0.5, after the outer sphere is accelerated instantaneously from fi 2 = 0.8 at t < 20 to 
f2 2 = 1 at t > 20, and the mutual friction is changed simultaneously from GM to HV. The 
snapshots correspond to (a) * = 21, (b) * = 22, (c) t = 50, (d) t = 100, (e) t = 120, and (f) 
t = 140. Time is expressed in units of Q^ 1 . 

Fig. 4. — Meridional streamlines of the normal fluid (left) and superfluid (right) components 
as a function of time for 5 = 0.5, AVI = 0.2, and GM mutual friction. There is no spin-up 
event at t — 20, unlike in Figure 3. The snapshots correspond to (a) t = 18, (b) t = 20, (c) 
t = 50, (d) t = 100, (e) * = 120, and (f) t = 140. Time is expressed in units of fif 1 . 

Fig. 5. — Fractional change in angular velocity of the outer sphere, AQ/Q = [fi 2 (t) — 
f2 2 (20)]/f2 2 (20), as a function of time, before and after a spin-up event at t — 20 where 
the mutual friction is changed instantaneously from GM to HV and Q 2 jumps according to 
f2 2 = 0.9 — > 1 (solid curve), f2 2 = 0.8 — > 1 (dashed curve), and f2 2 = 0.7 — > 1 (dotted-dashed 
curve). Time is measured in units of fif 1 . Dimensionless gap width: (a) 5 = 0.2, (b) 5 = 0.3, 
(c) 5 = 0.4, and (d) 5 = 0.5. 

Fig. 6. — Fractional change in the angular acceleration of the outer sphere, A(l/(l = [Q 2 (t) — 
Q 2 (20)]/J1 2 (20), as a function of time, before and after a spin- up event at t — 20 where the 
mutual friction is changed instantaneously from GM to HV and Q 2 jumps according to 
Q 2 — 0.9 — > 1 (solid curve), Q 2 = 0.8 — > 1 (dashed curve), and Q 2 = 0.7 — > 1 (dotted-dashed 
curve). Time is measured in units of Q^ 1 . Dimensionless gap width: (a) 5 = 0.2, (b) 5 = 0.3, 
(c) 5 = 0.4, and (d) 5 = 0.5. 
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Fig. 7. — (a) Evolution of the z component of the torque on the outer sphere (multiplied 
by 10 4 ) as a function of time for a superfluid with GM mutual friction (solid curve) and a 
classical Navier-Stokes fluid (dashed curve), with 6 = 0.5, Re = 3 x 10 4 , and Afi = 0.1. The 
torque is expressed in units of pR\Vt{ and the time in units of fl^ 1 . (b) Fractional change 
in the angular acceleration AQ/Q = [^(t) — ^2 (20)] /^2 (20) following the spin- up event 
tt 2 = 0.9 -> 1.0 at t = 20 with 5 = 0.5 and Re = 3 x 10 4 , in three cases: GM -> HV 
transition (dashed curve), superfluid with GM mutual friction (solid curve), and classical 
Navier-Stokes fluid (dashed-dotted curve). 



Fig. 8. — Meridional streamlines for superfluid spherical Couette flow with 5 = 0.5 and 
Re = 3 x 10 4 , obtained by integrating the in-plane velocity components in the plane x — 0. 
(a) Viscous fluid at t — 20, with fi 2 = 0.9 and Qi = 1.0, and (b) at t — 22, after a sudden 
spin up of the outer sphere Q 2 — 0.9 — > 1.0. (c) Viscous normal component at t — 20, with 
f^2 = 0.9, Q\ = 1.0, and GM friction, and (d) at t — 22, after a sudden spin-up of the outer 
sphere Q 2 = 0.9 — > 1.0, while simultaneously changing the friction from GM to HV. (e) 
Inviscid superfluid component at t — 20, with Q 2 = 0.9, Qi = 1.0, and GM friction, and (f) 
at t — 22, after a sudden spin up of the outer sphere f2 2 = 0.9 — > 1.0, while simultaneously 
changing the friction force from GM to HV. Time is measured in units of fi^ 1 . 



Fig. 9. — Turbulent-laminar vorticity transition during a glitch. Evolution of the torque on 
the outer sphere before and after the outer sphere is impulsively accelerated from f2 2 = 0.9 
to n 2 = 1-0 at t — 20, for 5 = 0.5. The angular velocity of the outer sphere is f2 2 = 0.9 at 
< t < 20, and Q 2 (t) = 1.0 + 0.001(t - 20) during the time interval 20 < t < 120. The five 
meridional slices are contour plots of \(v n )z/v-DG\ at (a) t = 20, (b) t = 22, (c) t = 50, (d) 
t = 70, and (e) t = 120. Dark regions indicate where the DGI is active, i.e. \(v n )z/vT>G\ ^ 1-0, 
white regions indicate where the DGI is not active, i.e. |(v n ) z /t>DG| < 1-0, and the mutual 
friction is of HV form. 
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